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1.0 Introduction 

The ability to accurately predict the thermomechanical deformation response of 
advanced composite materials continues to play an important role in the develop- 
ment of these strategic materials. Analytical models that predict the effective 
behavior of composites are used not only by engineers performing structural anal- 
ysis of large-scale composite components but also by material scientists in devel- 
oping new material systems. For an analytical model to fulfill these two distinct 
functions it must be based on a micromechanics approach which utilizes physi- 
cally based deformation and life constitutive models and allows one to generate 
the average (macro) response of a composite material given the properties of the 
individual constituents and their geometric arrangement. Only then can such a 
model be used by a material scientist to investigate the effect of different deforma- 
tion mechanisms on the overall response of the composite in order to identify the 
appropriate constituents for a given application. However, if a micromechanical 
model is to be used in a large-scale structural analysis it must be 1) computation- 
ally efficient, 2) able to generate accurate displacement and stress fields at both 
the macro and the micro level and 3) be compatible with the finite element 
method. Additionally, new advancements in processing and fabrication techniques 
now make it possible to engineer the architectures of these advanced composite 
systems. Full utilization of these emerging manufacturing capabilities require the 
development of a computationally efficient micromechanics analysis tool capable 
of accurately predicting the effect of microstructural details on the internal and 
macroscopic behavior of composites. The above mentioned computational effi- 
ciency is required since 1) the large number of parameters that must be varied in 
the course of engineering (or designing) composite materials, and 2) the optimiza- 
tion of a material’s microstructure will require the integration of the micromechan- 
ics model with optimization algorithms. From this perspective, analytical 
approaches that produce closed form expressions which describe the effect of a 
material’s internal architecture on the overall material behavior are preferable to 
numerical methods such as the finite element or finite difference schemes. 

A number of models presently exist that can fulfill some aspect of the aforemen- 
tioned tasks. However, there are very few working models that are both computa- 
tionally efficient and sufficiently accurate at the micro- as well as the macro-level. 
One such micromechanics model with the potential of fulfilling both tasks is the 
method of cells [1] and its generalization [2]. The comprehensive capabilities and 
efficiency of this method has been documented in references [4] and [5J. Conse- 
quently, the recently developed, computationally efficient and comprehensive 
micromechanics analysis code, MAC, who’s predictive capability rests entirely 
upon the fully analytical micromechanics model, herein referred to as the general- 
ized method of cells, GMC, [2 and 3] will now be described. MAC is a versatile 
form of research software that “drives" the double or triply periodic micromechan- 
ics constitutive models based upon GMC. GMC is capable of predicting the 
response of both continuous and discontinuous multi-phased composites with an 
arbitrary interned microstructure and reinforcement shape. GMC is a continuum 
based micromechanics model that provides closed-form expressions for the mac- 
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roscopic composite response in terms of the properties, size, shape, distribution, 
and response of the individual constituents or phases that make up the material. 
GMC also utilizes physically based viscoplastic deformation and life models for 
each constituent. Furthermore, expressions relating the internal stress and strain 
fields in the individual constituents in terms of the macroscopically applied 
stresses and strains are also available through strain or stress concentration fac- 
tors. These expressions make possible the investigation of failure processes at 
the microscopic level at each step of an applied load history. Similarly, GMC pro- 
vides the capability of studying the influence of bond strength at the fiber/matrix 
interface which recently has been shown to be an important damage mechanism. 

MAC enhances the basic capabilities of GMC by providing a modular frame- 
work wherein 1) various thermal, mechanical (stress or strain control) and thermo- 
mechanical load histories can be imposed, 2) different integration algorithms may 
be selected, 3) a variety of constituent constitutive models may be utilized and/or 
implemented and 4) a variety of fiber architectures may be easily accessed 
through their corresponding representative volume elements. Figure 1 illustrates 
the basic flow diagram for this modular framework. The capabilities for this version 
of MAC are discussed in section 2, whereas theoretical and background informa- 
tion on the basic capabilities itemized above is given in section 3. Section 4 
describes how one might use MAC and will be the most referred to section in the 
entire manual. Finally, section 5 gives some insights into the future modifications 
planned for MAC. 
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2.0 Current Capabilities 

In this section the current features/capabilities of MAC are itemized. 

• Load Types: 

- Thermal 

- Mechanical 

- Thermomechanical 

• Aboudi GMC Models: 

- Double Periodicity Model for continuous reinforcement 

- Triple Periodicity Model for discontinuous reinforcement 

• Graphical Output 

- Up to 5 x-y data plot files may be generated for both macro and micro (sub- 
cell) quantities 

- PATRAN for subcell geometry and color results evaluation, e.g. stress, strains, 
inelastic strains, J 2 , etc. 

• Integration Options: 

- Forward Euler 

- Predictor/Corrector 

• Constitutive Models 

- elastic orthotropic 

- inelastic viscoplastic models: 

- Bodner-Partom 

- Robinson 

- GVIPS 

- user defined model using subroutine USRMAT 

• Fiber/Matrix Interface Layer 

• Laminate Analysis 

• Fiber/Matrix Debonding Conditions 

• Free Format Input Data 
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3.0 Background 

3.1 Micromechanics Models 

As stated in the introduction, MACs predictive capabilities rest entirely upon the 
fully analytical micromechanics model known as GMC which is capable of predict- 
ing the inelastic response of both continuous (double periodicity) and discontinu- 
ous (triple periodicity) multi-phased composites with an arbitrary internal 
microstructure and reinforcement shape. Prior to describing the available architec- 
tures (Representative Volume Elements, RVE’s) within MAC as discussed in sec- 
tion 4, a brief overview of the theoretical foundation behind the generalized 
method of cells follows. A more complete discussion of the theoretical formulation 
is given in [1] - [3]. 

In the original formulation of the method of cells, a continuously (or discontinu- 
ously) reinforced, unidirectional fibrous composite is modeled as a rectangular, 
double-periodic (or triply-periodic) array of fibers embedded in a matrix phase. 
The periodic character of the assemblage allows one to identify a repeating unit 
cell that can be used as a building block to construct the entire composite. The 
properties of the repeating cell are thus representative of the properties of the 
entire assemblage. The unit cell consists of a single fiber subcell surrounded by 
three matrix subcells for continuous and seven for discontinuous composites, 
hence the name method of cells. The rectangular geometry of the repeating unit 
cell allows one to obtain an approximate solution for the stresses and strains in 
the individual subcells given some macroscopically homogeneous state of strain 
or stress applied to the composite. The approximate solution to the posed bound- 
ary value problem is, in turn, used to determine macroscopic (average) or effec- 
tive properties of the composite and the effective stress-strain response in the 
inelastic region. 

In the generalized method of cells for continuous (or discontinuous) fibrous 
composites, the repeating unit cell can consist of an arbitrary number of phases. 
Hence the generalized method of cells is capable of modeling a multiphase com- 
posite. This generalization extends the modeling capability of the original method 
of cells to include the following: 1) inelastic thermomechanical response of mul- 
tiphased metal matrix composite, 2) modeling of various fiber architectures 
(including both shape and packing arrangements), 3) modeling of porosities and 
damage, and 4) the modeling of interfacial regions around inclusions including 
interfacial degradation. 

The basic homogenization approach taken in the micromechanical analysis 
consists essentially of four steps. First, the repeating volume element, RVE, of the 
periodic composite is identified. Second, the macroscopic or average stress and 
strain state in terms of the individual microscopic (subcell) stress and strain states 
is defined. Third, the continuity of tractions and displacements are imposed at the 
boundaries between the constituents. These three steps, in conjunction with 
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micro-equilibrium, establish the relationship between micro (subcell) total, thermal 
and inelastic strains and macro (composite) strains via the relevant concentration 
tensors. In the fourth and final step, the overall macro constitutive equations of the 
composite are determined. These four steps form the basis of the micro-to-macro- 
mechanics analysis which describe the behavior of heterogeneous media. The 
resulting micromechanical analysis establishes the overall (macro) behavior of the 
multi-phase composite and is expressed as a constitutive relation between the 
average stress, strain, thermal, and inelastic strains, in conjunction with the effec- 
tive elastic stiffness tensor. 

That is, 


g = B*(e-£ -f) CEQ1) 

where for the most general case of discontinuous reinforcement with N a by N p by 
N y number of subcells, the effective elastic stiffness tensor, B , of the composite 
is given by, 


N a N p N 7 

B - =miLLL Wt? 4 

a=lp = ly=l 

the composite inelastic strain tensor is defined as, 


-i 

£ = 


-B 


*- 1 


N„ 


N* N y 


dhl 


E E E 

a - ip = ly = l 




c (aPY) ( ^«*PY) g /_|/(«PY) ) 


the average thermal strain tensor as, 


(EQ2) 


CEQ3) 


-T 

£ 


a= ip = ly = 1 


(EQ4) 


and e is the uniform applied macro (composite) strain. For the case of continuous 
reinforcements with N p by N y number of subcells, eq. (2) - (4) reduce to the follow- 
ing: 


N, N y 


B ~* = hi* I VtC (Py V P y) 

(EQ5) 

P=lY=l 



(EQ6) 

P- it- i 


*_1 N, N y 

B r (PY) m (P y)T -r(PYh 

hi I LVy^ §s-§ ) 

(HQ 7) 


P = 1 Y = 1 
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In the above equations matrix notation is employed; where, for example, the 
average stress, g, average applied strain, f , and inelastic subcell strain, , vec- 
tors represent, 


g - { a n> °22’ a 33> a 12» a 23* CT 13^ (EQ 8) 

I = {^ll’^22’^33> ^12>^23’^13^ (EQ9) 




_/(/v a /v„/v T ) 

•J § J 


(EQ 10) 


where the six components of the vector e /(aPY) are arranged as in eq. (9). Similar 
definitions for §' , | /(aPY) also exist. Note that the key ingredient in the construc- 
tion of this macro constitutive law is the derivation of the appropriate concentration 
matrices, A ( “ Py) and D ( “ Py) having the dimensions 6 by 6 and 6 by 
respectively, at the micro (subcell) level. The definitions of A and D, although not 
given here, may be found in references [2] and [5]. Finally, the matrix C a Y rep- 
resents the elastic stiffness tensor of each subcell (aPy) and d a , /i p , l y the 
respective subcell dimensions (see Fig. 2) wherein, 


N„ 


N* 


N, 


d = 


= E 


h = 


a = 1 


E' 

p = l 


/ = 


■ E 


y=l 


Similarly, given the concentration matrices A (<xPy) and £> ( “ Py) , expressions for 
the average strain in each subcell can be constructed, i. e., 


i (oPT) = A (apY) | + D (aPY) (e^ + e^) 


as well as average stress, 


^(a P Y) = g , ( aPY ^ [i 4 (aPY) | + ^^ aPY) (g^ + §j ) - (e^ a ^ Y) + § r(a ^)] 


The analytic constitutive law, see eq. 1 , may be readily applied to investigate the 
behavior of various types of composites, given knowledge of the behavior of the 
individual phases. Numerous advantages can be stated regarding the current 
macro/micro constitutive laws as compared to the other numerical micromechani- 
cal approaches in the literature, e.g. the finite element unit cell approach. One 
advantage is that any type of simple or combined loading (multiaxia! state of 
stress) can be applied irrespective of whether symmetry exists or not, as well as 
without resorting to different boundary condition application strategies as in the 
case of the finite element unit cell procedure. Another, advantage concerns the 
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Double 

Periodicity 




Figure 2: Subcell Dimension Nomenclature 
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availability of an analytical expression representing the macro elastic-thermo- 
inelastic constitutive law thus ensuring a reduction in memory requirements when 
implementing this formulation into a structural finite element analysis code. Fur- 
thermore, this formulation has been shown to predict accurate macro behavior 
given only a few subcells, within the repeating cell (see references [2], and [4]). 
Whereas, if one employs the finite element unit cell procedure, a significant num- 
ber of finite elements are required within a given repeating unit cell to obtain the 
same level of accuracy as with the present formulation. Consequently, it is possi- 
ble to utilize this formulation to efficiently analyze metal matrix composite struc- 
tures subjected to complex thermomechanical load histories. This is particularly 
important when analyzing realistic structural components, since different loading 
conditions exist throughout the structure, thus necessitating the application of the 
macromechanical equations repeatedly at these locations. 

3.2 Integration Algorithms 

There are two integration algorithms currently available within MAC. The first is 
the standard, explicit Forward Euler algorithm, which can be expressed as, 

w' 1 + 1 = w, + Ar/(r l , w,.) 


where 


= y(t t ) 

w «+ i = y('« + i) 


and f{t it w,) is the rate of change with respect to time of the vector, y(r f ) , i.e. 


/('«> w,-) 


= y = 


dyU,) 

dt 


The second is a predictor/corrector algorithm which uses: 
1) a 4 th order Runge-Kutta starter: 


w i = w ,-i + 


(k 1 + 2K 2 + 2k 2 + k a ) 


where 
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ATi = A f/Xfj.i.w,--!) 

K i = A'/(<i - 1 + Y» l *',-i + y) 

*3 = A //((,. 1 + y , H ',- i + y ) 
^4 = Atf(tj_ 1 + At, Wj_ j + Kj) 


2) with an Adams Bashforth four-step predictor: 




p 

W ■ i = 
I + 1 


= a 1 Wj = a 2 h> 2 = oc 3 w 3 = a 4 
+ ^[55/(r„ w i )-59f(t i _ v w i _ 1 )+ 37/(f, _ 2 , _ 2 ) 


-9 /(?,-_ 3 . W i-3>1 


and 3) an Adams Moulton four step corrector: 


W; 


Ar 


f + i = w, + g[9/(f l + i.< + i) + l 9 /(f l .w i )-5/(r i _i,w i . 1 ) 


+ /(**- 2* w i-2)l 


where the a’s come from the 4 th order Runge-Kutta starter. Further details may 
be found in [6]. 

It has been found, based on experience, that for relatively rapid monotonic or 
cyclic loadings it may be more efficient to use the Forward Euler integrator since 
the predictor/corrector requires 5 evaluations per step, as shown above. However, 
in the case of creep, relaxation or slow monotonic or cyclic loading histories, sig- 
nificant increases in solution speeds can be obtained using the predictor/corrector 
algorithm with a self-adaptive time step. 
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Finally, within MAC the vector y as used above contains the following macro 
quantities: 


position 

contents 

1 -6 

Macro Total Strain 

7 - 12 

Macro Stress 

13-18 

Macro Inelastic Strain 

19-30 

(currently empty 
space for 2 6x1 vectors) 

31 -36 

Macro Thermal Strain 

37 

Current Temperature 


and the quantities associated with each subcell are stored sequentially in y, such 
that 


position 

contents 

38-43 

Micro Total Strain 

44-49 

Micro Stress 

50-55 

Micro Inelastic Strain 

56-67 

Micro Internal State Variables 
(space for 2 6x1 vectors) 

68-73 

Micro Thermal Strain 


The above 36 positions are repeated for the total number of (N) subcells thus 
bringing the total length of the y vector to 37+36N. It follows that a second vector 
of similar length contains the corresponding macro and micro rates y. 

3.3 Available Constituent Constitutive Models 

Currently MAC provides one elastic and four inelastic constitutive models. 
These models have been selected purely based upon the availability of material 
parameters for the materials of interest. However, MAC is designed in a modular 
fashion thus allowing the implementation of additional constitutive models through 
a user defined subroutine. Two of the four available inelastic models are capable 
of representing transversely isotropic material behavior, thus allowing one to 
investigate the reinforcement of an anisotropic matrix or the idealization of an 
anisotropic fiber. In all four inelastic models a purely elastic response is possible 
by modifying a single material parameter for each model as noted below. 
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3.3.1 Elastic Model 

Reference: Mechanics of Composite Materials, Jacob Aboudi, Elsevier, 1991 

The following transversely isotropic model is provided for those materials that 
have an elastic only response, such as the fiber constituent in a composite. 


G 11 


Cu 

C 12 

C 12 

0 

0 

0 


e ll 

a 22 


C 12 

C 22 

G 23 

0 

0 

0 


e 22 

a 33 


C \2 

C 23 

C 22 

0 

0 

0 


633 

g 23 


0 

0 

0 

C44 

0 

0 


Y 23 

ct 13 


0 

0 

0 

0 

C44 

0 


T 13 

_ a 12 


0 

0 

0 

0 

0 

C G6 


Y 12 


where components C tJ can be expressed in terms of five independent constants, 


thus, 


E a , E t , w a , V r , G a 


C n = E a + 4kv a 
C l2 = 2kv a 

0.5E t 

C 22 = K+ (TT^ 


C')'! — K- 


0.5Ei 


23 " (l+v r ) 

g a 

( G 22 ~ G 23 ^ 


C A4 = g a 


C 66 ~ 
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3.3.2 Bodner-Partom Model 

Reference: Mechanics of Composite Materials, Jacob Aboudi, Elsevier, 1 991 

This model represents the Bodner-Partom viscoplastic model with isotropic 
hardening, Z, and can be used for an initially isotropic metallic material. 

The flow law is given as: 


where 



i 



n 



S ij ~ ®ij 


The evolution law for isotropic hardening is given as: 

W PL 

Z = m(Zi-Z CJ y)— 

. PL . 

where Z 0 , Z Y and m are inelastic constants and the plastic work rate, W , is 
given by; 



°ifij 


t 3 t 

Z eff = Z 0 + <?|Z (t)^t + (1 -q) £ r j 7 JZ(T)r jy (T)dx 

o «. i = i o 

r,j(0 = a, v (t)/[a„(I)<J H «)] 1/2 


An elastic only response may be obtained by setting the material parameter D 0 to 
zero. 
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3.3.3 Robinson Creep Model 


Reference: NASATM 103172, 1990 

This model represents a transversely isotropic material, wherein the vector of 
direction cosines, d { defines the strong material direction along which no inelastic- 
ity occurs. 

Flow Law: 


./ 



£o 



CT o 


Jf 


Evolution Law: 


Cl; 


— RG m ~ p ^ 
G P2 ^o OqJg 


where 

r ij = Ly-^KiDy-bu) 

n,y = a v -\H3D,j-b v ) 


in which 


F = 4 * 

c o 

♦ = J 2-\* 2 

4 

v = h--f 

ICS o 

1 b 

n 


and 


J 1 — 




.. _ i 

— 2 a ij a ji 


D^ji &ij ~ ^idj ^ ij Sij a ij 

I = Dydji Sjj = Ojj - 


Special cases involving an isotropic material can be obtained by defining 
d i = (1/3, 1/3, 1/3) and/or elastic only response by setting a 0 to an extremely 
large number. 
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3.3.4 Robinson Viscoplastic Model 
Reference: NASATM 105819, 1992 

This model represents a transversely isotropic material wherein the vector of 
direction cosines d i defines the preferred material direction. In this model the 
strength of anisotropy is specified by the parameters to and rj ; where co is the 
ratio of the normal longitudinal and transverse yield stress and r\ is the ratio of 
longitudinal and transverse shear strengths. 


Flow Law: 


Evolution Law 


./ 

e U 


_ <F n ) 


2\l ^ 




n. 


where 


and 


r y - ^ij - K^ki^jk + Djk^-ki 2^ /() ^ 3£> V 

n ij = a (j - fyDfodji + D jk a kl - 21 - ^1 q(3D-j - 5,y) 


= \\l x + — 2 i 2 + \ / 3 1 - 

T| 4(4co - 1) J 

G = m + l/2 + ^ ^1 

4 L r\ 4(4co - 1) J 

G = (G - GQ)Hv[S i jK i j] + G 0 


h = 


1 

II 

CM 

^2 ~ 2 

1 = 

D„ = d i d j 



r _ 4(co 2 - 1) 
*» 2 


n 

4 m - 1 


16 



MAC M kromecharecs A nalysis C ode 


The invariants 1 1 , / 2 , 7 3 are the same as those given above but with L- replaced 
by a,». Special cases involving an isotropic material and/or elastic only response 
can lae obtained by defining to = r) = 1 and/or by setting k t to an extremely large 
number. 
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3.3.5 Generalized Viscoplastic Potential Structure (GVIPS) Model 
Reference: NASATM 106609, 1994 

This model is a fully associative, multiaxial, isothermal, nonlinear kinematic 
hardening viscoplastic model for use with initially isotropic metallic materials. A 
unique aspect of this model is the inclusion of non-linear hardening through the 
use of a compliance operator Q^ kl in the evolution law for the back stress. This 
non-linear tensorial operator is significant in that it allows both the flow and evolu- 
tionary laws to be fully associative (and therefore easily integrated) and greatly 
influences the multiaxial response under non-proportional loading paths. 

Flow Law: 


J = SL'lik 

“ 2 


if FZ 0 


where 



Internal constitutive rate equation 

@ij ~ Lij rs A rs 


Evolution Law: 


An = 4 - - iRaB f “n if 


~ Qrslm E lmnp np 


3R a B 0 G q 


•a 


np 


a ij E ‘ij ^ 0 


if afa<0 


where 
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F = {&-Y) 

K 

Y = < 1 - P7g> 



K 


'yrs 


= IQijrsV = 


K 


35 0 (l+5 1 pG p ' 1 )'. 


^ijrs ~ 


3S 1 (p-l)G / 


K 2 (l + fl 1 pG p_1 (6/»-5)) 


a rs a ij 


and 

h = f a y fl y L ij = S (T a i/ 

J 2 - 2 ^y^y ^y = °ij~ 2°klfiij 

The special case of an elastic only response maybe obtained by setting k to an 
extremely large value. 
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3.4 Laminate Theory 

Version 2.0 of MAC now includes the capability to analyze composite laminates. 
In-plane axial displacements and membrane forces may be applied to the lami- 
nate. That is, the global laminate stress-strain relation that is solved within MAC is 
expressed as, 


p 


— 


- 


- - 




—7 

Nxx 


^11 

A 12 

A 13 


^XX 


Nxx 


Nxx 

Nyy 

= 

A 21 

A 22 

^23 


£yy 

- 

n'yy 


N T yy 

N X y_ 


a 31 

A 32 

a 33 


5*2 


F'xr 


n t xy _ 


y = A£-N ! -N T (EQ 12) 

where N , N l , and N T are the global laminate total, inelastic and thermal stress 
resultants, respectively. The matrix A is the global laminate stiffness and § the 
global laminate strains. 

In forming the laminate stiffness A the generalized method of cells model, 
GMC, is utilized to calculate the individual lamina properties. In this regard, the 
individual laminate stiffness, in lamina coordinates, Q is given by, 



Qn Qn 0 
021 0-22 0 
0 0 <266 


in which the components of Q are given as, 


Qn 

Q 22 


C 11- 
C 22 “ 


^13^31 

C 33 

^23^32 

C 33 


Ql2 ~ C \2 ~ 
Q33 = C 


C\-iC 


13 w 23 


'33 


66 


(EQ 13) 


(EQ 14) 


The C i} in the above are the effective macro properties for the unidirectional com- 
posite lamina obtained from GMC. 

It follows that the lamina stress-strain relation in global (laminate) coordinates is 
given by the relation, 
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or 


where 



Qn Q12 

Q 12 022 

Q 13 @23 


Qn 

<223 

Q33 


kf x Z 


g = Qe 


Q = R l QR 


(EQ 15) 


(EQ 16) 


(EQ 17) 


R = 


COS0 2 

sin0 2 

-sin0cos0 


2 

sin0 2sin0cos0 

cos0 2 (-2)sin0cos0 

2 . 2 
sin0cos0 cos0 - sin0 


(EQ 18) 


and 0 is the orientation of the longitudinal lamina axis with respect to the global x 
direction, thus Q k is the transformed lamina stiffness, i.e. from local lamina to 
global laminate coordinates. In addition, g and | are the lamina stress and strain 
in laminate coordinates. It then follows that the global laminate stiffness A is given 
by, 


nly 

4=1 2A ® 19 > 

*= 1 

in which nly is the total number of layers in the laminate and t k is the thickness of 
the k th lamina. 

Returning to EQ. 12, the quantities N 1 and N T (the laminate inelastic and ther- 
mal stress resultants, respectively) are calculated from the individual lamina con- 
tributions through the following relations, 

nly nly 

f = £jy/r t jv r = £ N k T t t (EQ20) 

k= 1 * = 1 

I y 

The individual lamina stress resultants N k and N k are given by, 

N k = Q k e k (EQ21) 
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3.5 User Defined Material Model 

Another addition to version 2.0 of MAC , is the option for a user to implement 
their own constitutive model. This is accomplished by using the subroutine 
USRMAT into which the user writes the necessary FORTRAN code for the partic- 
ular constitutive model being implemented. The USRMAT subroutine is as shown 
below. 


USRMAT 



SUBROUTINE USRMAT(DSA f SA,PE,PV, 

& TIME,TSTERCTEMRDTEMPR,TREF,NIO,NPE,NPV,NCE) 


c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


purpose: user supplied constitutive model 

note: 1) in this subroutine, [sa] and [dsa] contain the 

“micro” quantities for aboudi’s micromechanics model 

2) arrangement of [dsa] & [sa] arrays are as follows 


variable 

location 

1 

| strain rate 

(1-6) 

j stress rate 

(7-12) 

j inelastic 
j strain rate 

(13-18) 

1 12 “slots” 

(19-30) 

| for state variables 

|— 

| thermal strain rate (31 -36) 

1 


IMPLICIT DOUBLE PRECISION (A-H.O-Z) 
c 

DIMENSION DSA(36),SA(36),PE(NPE),PV(NP V) 

c start user supplied material routine c 

c 

c place code defining model here 
c 

c end of user supplied material routine c 


RETURN 

END 
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A description of the input and output required for the usrmat subroutine is as fol- 
lows 

Program data available to USRMAT: 

SA array containing current total quantities for all 

of the state variables 

PE array containing elastic constants 

PV array containing inelastic constants 

TIME current time 

TSTEP current time increment (step) 

CTEMP current temperature 

DTEMPR current temperature rate 

TREF reference temperature 

NPE total number of elastic constants 

NPV total number of inelastic constants 

NIO output file unit number 

NCE current subcell number 

Output expected from USRMAT: 

DSA current increments in all state variables 

The state variables are arranged in SA and DSA in the following order: 

Position Quantity 

1 -6 strain 

7-12 stress 

13-18 inelastic strain 

1 9-30 available space for 2 6x1 vectors for model 
dependent internal state variables 
31-36 thermal strain 

Again, SA contains the total quantities and DSA contains the rates- 

is- Note: Appendix H contains a sample input file and the USRMAT subroutine 
containing an implementation of the Bodner-Partom viscoplastic material 
model. 
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3.6 Interface Modeling 

Interfaces in composite materials play a major role in the determination of their 
mechanical and thermal properties. Consequently, it is important to have the abil- 
ity to model interface behavior accurately. This is accomplished in MAC/GMC in 
one of two ways. The first is to define an actual interface region with its own con- 
stitutive behavior. In this way the influence of initial imperfections (flaws, voids, 
improper wetting, etc.) and induced interfacial damage (due to stress, environ- 
ment, chemical reactions, etc) may be incorporated into the micromechanical 
analysis of the overall behavior of the composite. The development of proper inter- 
facial constitutive models is an active area of research and MAC , through the use 
of its USRMAT routine, provides the researcher with a convient tool for testing new 
and existing interfacial constitutive models. 

The second approach to modeling the effect of imperfect (weak) bonding 
between two phases (e.g. a fiber and matrix) is to assume a jump in the displace- 
ment field at a particle-matrix interface may occur given certain conditions, while 
still maintaining continuity of the traction vector. In the spirit of Jones and Whittier 
(see ref .10) and Achenbach and Zhu (see ref. 1 1) we have assumed the following 
“elastic-perfectly plastic” flexible interface model. 


(uj = R n • aj) 
(u/ = R, ■ a/) 


if 


°n ^ °DBn 
I . 

°t ^ °DBt 


(EQ 23) 


<“/ = o) 
(«/ = 0 ) 


a n < °DBn 

a/ < o DBl 


(EQ 24) 


where R n M t , o DBn and o DBt are the interfacial normal and shear , strength and 
debond stress, respectively. Implementation of this model impacts the definition of 
the concentration matrices, A discussed in section 3.1. 
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4.0 Running MAC 

MAC has been developed for a Sun IPX workstation running Sun Solaris. Upon 
linking the associated object modules that comprise MAC, execution is begun by 
simply typing MAC (in upper case) and hitting return. 

4.1 Input and Output Files 

MAC requires that the names of the input and output files be specified. The user is 
prompted for the names of the input, infilename, and output file, outfilename, and 
the names of the input and output files maybe anything the user desires. As will be 
explained in the following sections, the user will also need to specify the names of 
the PATRAN output files and the files containing the X-Y data for plotting. 

Note that a “debug” file may be generated depending on the “PRINT” level the 
user specifies in the input file (see 4.2.2). This file will have the same name as the 
output file but with the string “_debug” automatically appended at the end of the 
filename, i.e., outfilename_ debug. 

Input Output 

infilename outfilename 

debug file ( outfilenamejiebug) 

PATRAN files (optional, see section 4.2.11) 
plot data files (see section 4.2.15) 

4.2 Input Requirements 

This section describes how the data should appear in the input file, infilename. 
In the following, each block of input data will have its own subsection and will typi- 
cally contain the following information: 

1) statement of purpose 

2) declaration of input data block 

3) example(s) 

4) notes 

A data block has the following general format: 

1 ) ^keyword beginning of data block 

2) input line(s) 

3 ) % end of data block 
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where 

1 ) * keyword: (denotes beginning of data block) 

Each input data block begins with the appropriate keyword, which starts with 
the * symbol. For example, the keyword for selecting the method of integration 
is *SOLVER. The input routine scans the input file and locates the appropriate 
keyword and then reads the corresponding input data 

2) Input line(s): 

The input lines contain the necessary input data. The specific format of these 
input lines will be given in the next section. 

Two special characters (& and #) are provided for entering the input. 

The (continuation) symbol: 

For input data that is too long to fit on a single line, the symbol is used 
as the continuation character. Thus, a line of data may be divided into a 
series of lines. For example, when entering the material properties all of 
the data cannot fit on a single line, thus the continuation character is 
required: 


EL=1 1 700., 1 1 700. ,0.365,0.365,4287.5,1 .,1 . & 

VI=0.8E-8,0.1 ,0.1 E-5.0..0.85E-3, 0.05,1 .,1 .,1 .,3.3,1 .8,1 .35,1 .,0.01 

The “#" (comment) symbol: 

The input file may also contain “comments" for the users convenience. The 
“#” symbol is used to mark a comment line, with the requirement that the 
“#” symbol appears in the first column. 

3) % (denotes end of data block) 

Each input block ends with the “%” symbol. The % symbol must b£ included 
as it signals the input routine that a particular data block has been completed. 

A sample data block is as follows, 

*SOLVER 

# NTF=2 ISTM=0.0001 ERR=0.1E-3 
NTF=2 ISTM=0. 0000024 ERR=0.1E-2 % 

Note how the # symbol is used to comment out a line of input data. Thus, the 
user can change various parameters by simply commenting out input lines con- 
taining different input data. Also note how the data block is terminated with the % 
symbol and how it need not appear on a separate line. 
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es* Note: It is suggested that the order of keywords in a given input file follow that 
given in this manual.: 


4.2.1 Header Line: 

Purpose: Define the title of this particular job (80A format). 

Note: the Header Line is limited to one 80 character line, and is always taken 
as the first line in the input file. : 

problem title 

Example: Transverse tensile response of 35% SCS6/Ti--6-4 


4.2.2 Output Print Level : 

Purpose: To control the output generated 

*PRINT 

NPL=np/v/ 


Where 

nplvl: 

-1 = print out macro stiffness matrix, engineering constants and macro ther- 
mal expansion coefficients 
0 = minimal output (most commonly used) 

7 = basic program trace 

10 = program execution trace and all array data 

WARNING: this generates a very large output file 

Example: (minimal print out) 

*PRINT 

NPL=0 % 
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4.2.3 Load Type : 

Purpose: To select load type 

*LOAD 

LCON=nse/ LO P=lop LSS =iopt 


% 

Where 
nse I: 

1 = Thermal Load 

2 = Mechanical Load 

3 = Thermomechanical Load 

lop: (for 2-D and 3-D GMC) 

1 = axial load in 1 -direction 4 = shear load 23-direction 

2 = axial load in 2-direction 5 = shear load 13-direction 

3 = axial load in 3-direction 6 = shear load 12-direction 

For 2-D GMC : 

7 = axial loads in 1 and 2 directions 

8 = axial loads in 2 and 3 directions 

9 = axial load in 1 -direction and shear load in 23-direction 
10 = axial load in 2-direction and shear load in 13-direction 

For 3-D GMC : 

7 = axial loads in 1 and 2 directions 

8 = axial loads in 2 and 3 directions 

9 = axial loads in 1 and 3 directions 

10 = axial load in 1 -direction and shear load in 23-direction 
11= axial load in 2-direction and shear load in 13-direction 
12 = axial load in 3-direction and shear load in 12-direction 


iopt. 

1 = Strain control 

2 = Stress control 
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4.2.4 Mechanical Load Control: 

er Note: This block Is only required if LCON = 2 or 3 

Purpose: Select type of load control for mechanical load 

*MECH 

NPTW=npfw TI=f 7 ,f 2 -.-fnpfw LO = h^-^nptw 


% 

Where 

nptw: number of points on load curve 
tp ■•■tnptw • values 
Ip I 2 , —I nptw •' ,oad curve values 

Example: 

3HVIECH 

NPTW =3 Tl=0., 1., 2. LO=0., 20., 30. % 



«*■ Note: For LOP=7,8,...etc (i.e. two load/displacement components) the following 
format is required: 

*MECH 

NPTW =2 Tl=0., 1.5 LO=0. P 0.015 component 1 curve 

NPTW=2 Tl=0., 1.5 LO=0., 0.01 % component 2 curve 


(See Appendix G for an example) 
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4.2.5 Temperature Control: 

«*■ Note: This block is only required if LCON = 1 or 3 

Purpose: Select control for temperature 

♦THERM 

NP l I =nptt Tl = tpt 2 — tpptw TE= — t&nptt 

% 

Where 

nptt: number of points on temperature curve 

t P t 2 , ..t nptw : time values 

top t& 2 , -t&nptt '■ temperature curve values 

Example: 

♦THERM 

NPTT=3 Tl=0., 1., 2. TE=100., 200., 150. % 
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w Note: Forth© Thermomechanical Load t np(w mt nptt , and both curves must have 
= o. But the number and time value of the data points in-between maybe dif- 
ferent, see figure 4. 



Figure 4: Load History Specification 
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4.2.6 Reference Temperature: 

Purpose: Define initial (starting) temperature 

*TREF 

TREF=fre/ 


% 

Where 

tref: 

— temperature 

Example: (set reference temperature at 70.) 

*TREF 

TREF=70. % 

*a* Note: only important in thermal or thermomechanical load case when using 
secant values of coefficient of thermal expansion. 

ns 3 Note: Units must be consistent with selected material model parameter units 
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4.2.7 Micromechanlcs Model Identification: 

Purpose: Select desired GMC micromechanics model 

*MODEL 

MOD =modid MATSYS=mateys 

ra- Note: If modid=3 (laminate option) you must enter the following line of data, 
NLY=n/y THK=thk 1t thk 2 ,... thk n , y CON =c 7> c 2 ,...c niy SYS =s 1 ,s 2 ,...s nly 
ANG=a 7 , 32,... 4n/y 

% 

Where 

modid: 

1 = double periodicity 

2 = triple periodicity 

3 = laminate option 

matsys: - number of different material systems in the analysis 
Note: matsys is required only for lam inate option (MOD=3) 

nly. - number of layers in the laminate 

thkf, thk 2 ,... thk nly : - thickness of each layer 

Ci, c 2 , ... c n i y : - layer material type id: 1 = isotropic 2 = anisotropic 

s 7 , s 2 , ... s n i y : - material system id number 

a 7 , 32, ... a n /y : - angles of each layer 


Example 1 : 

double periodicity model 

*MODEL 

MOD=1 % 

Example 2: 

laminate option with 1 material system 

JfrMODEL 

MOD=3 MATSYS=1 NLY=4 THK=0.25, 0.25, 0.25, 0.25 CON=2,2,2 I 2 & 
SYS-1 ,1,1,1 ANG=45.,-45.,-45.,45. % 
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4.2.8 PATRAN Output: (Optional) 

Purpose: To select PATRAN output file generation. 

*PATRAN 

FN ^prefix STP= npstp STO =Si,s 2 ,.-s nps p 


Where 

prefix'. 

— this will be the filename prefix that will be assigned to all PATRAN files. 
npstp: 

the number of steps to be stored for plotting 

Si,S2,--S n p S p- 

step numbers at which data is to be stored 

is* Note: This is an optional input line. If no PATRAN files are desired, do not 
include this data block. N 

Example 1 : (select PATRAN output, use file prefix “runl”) 

*PATRAN 

fn=run1 STP=3 STO=1 ,50,100 % 

«ar Note: If this option is chosen, MAC generates the following PATRAN files: 
(all files are written in PATRAN 2.5 neutral file format in which the subcells 
are treated as “elements”) 

1) prefix. patgeo: contains “geometry” information of the RVE. 

2) prefix. patstr: contains micro-stress quantities 

(file is formatted as a PATRAN 2.5 element results file) 

3) pref/x.patepsin: contains micro-inelastic strain quantities 

4) pre f/x.patepsto: contains micro-total strain quantities 

For Example 1 , the following files would be produced: 

runl .patgeo 
runl .patstr 
run1.patepsin 
runl.patepsto 

and will contain data for load steps 1 , 50, and 1 00 
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4.2.9 integrator Identification 

Purpose: Select type of integration scheme 

*SOLVER 

NTF=ntf 

If NTF=1 , enter the following line 

NPTS=nptS TIM= ■■■^npts STP=Sfj, St£, • • - Scripts- 1 

If NTF=2, enter the following line 

I STM =/s/pm ISTT =ispt ERR =errtol 


% 

Where 

ntf: 

1 = Forward Euler method 

2 = Predictor/Corrector method 

istpm : 

= initial mechanical load time step 

istpt: 

= initial thermal load time step 

» Note: only istpm is required for a mechanical load (nseh 2) , only istpt is 
required for a thermal load {nseh 1) and both istpm and istpt are required 
for a thermomechanical load {nseh 3) 


Stj 

= initial time step for that load increment 
errtol: 

= error tolerance for predictor/corrector 


«y Note: errtol is only required when using predictor/corrector (ntf = 2) 
Suggested errtol = 0.1 


Example 1: for a mechanical load, select predictor/corrector, with initial time 
step = 0.001 and error tolerance = 0.01 


*SOLVER 

NTF=2 

ISTM=0.001 ERR=0.01% 
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Example 2: select forward euler, time step = 0.001 

*SOLVER 

NTF=1 

NPTS=3 TIM=0.,1.,2. STP=0.01,0.02 % 


sf 7 = 0.01 


st 2 = 0 . 02. 


L 

ti= o. 



^Tlme 

t 3 = ' 2 - 


4.2.10 Thermal Conductivity Calculation: (Optional) 

Purpose: Calculate thermal conductivity of the composite 

*COND 

is? Note: If the *COND keyword is not used, the thermal conductivity calculations 
are ignored. 
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4.2.11 Constituent Material Model Identification: 

Purpose: To select the model for the fiber and matrix constituents 

*FIBER 

NFIBS=n//bs 


The following line is to be repeated for each fiber (nfibs): 

NF=n// MS =msf MF =ncmdf NDPT =dpt 

\DB= idbf f =1,2,... nfibs 


^MATRIX 

NMATXsnmatx 


The following line is to be repeated for each matrix (nmatx): 

NM=nm m MS =ms m MA =ncmd m NDPT=dpf MPJ=mat m 
IDB=/db m m=1,2,...nmatx 

*MONOL (available for laminate option only) 

blMON-nmon 

The following line is to be repeated for each isotropic material (nmon): 

NMO =nmo iso MS =ms iS0 MMO =ncmd iso NDPT=ojof MfiJ-mat iS0 
IDB =idb iso iso=1,2,...nmon 


% 

Where 

nfibs: 

number of different fibers 
nmatx: 

number of different matrices 
nmon: 

number of different isotropic layers 


nff. 

fiber material designation number running from 1 to nfibs, sequentially. 


nm w : 

matrix material designation number running from 1 to nmatx, sequentially. 
nmOj S0 . 

isotropic material designation number running from 1 to nmon. 


ms/. 

fiber material system id (required only when using laminate option) 
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ms m- . v 

matrix material system id. (required only when using laminate option) 


ms iso ■ ... 

isotropic material system id (required only when using laminate option) 

ncmd: 

material model identifier for either fiber, matrix or isotropic layer: 

1 = Bodner-Partom Model 

2 = Robinson Creep Model 

3 = Robinson Non-normalized Viscoplastic Model 

4 = Generalized Viscoplastic Potential Structure (GVIPS) Model 
6 = Elastic Material Model 

99= User defined model (see note on page 44 for special format instructions) 
mat: 

material identification letter for either fiber or matrix, selected from material 
database, see Table I or II depending upon. 

idb: 

read material constants from the database (y) or read from input file (n) 

dpt (This flag is only required for the Bodner-Partom model) 

flag indicating whether material constants should be temperature independent 
or temperature dependent 

1 = Temperature Independent 

2 = Temperature Dependent 

Table I: Isothermal Material Constants 

MAC Fiber Database 


Model 

Material 

Units 

mat 

Linear Elastic 

Boron 

MPa 

A 

ncmd = 6 

SCS-6 

MPa 

B 


Tungsten 

ksi 

C 


SCS-6 

ksi 

D 
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MAC Matrix Database 


Model 

Material 

Units 

mat 

Bodner-Partom 
ncmd = 1 

Aluminum 

MPa, sec 

A 

Robinson Creep 
ncmd = 2 

Kanthal 

ksi, hr 

A 

Robinson Viscoplastic 

Kanthal 

ksi, hr 

A 

ncmd = 3 

FeCrAlY 

ksi, hr 

B 

GVIPS 
ncmd = A 

TIMETAL21S 

ksi, sec 

A 


Example 1 : 

select 1 fiber, Boron, and 1 matrix material, Timetal21S, both read 
from database 

the Boron Timetal21 S composite 


♦FIBER 

NFIB=1 

NF=1 MF=6 MAT=A IDB=y% 

♦MATRIX 

NMATX=1 

NM=1 MA=4 MAT=A IDB=y% 

ar Note: See Section 3.3 for a mathematical description of each material model. 

Example 2: (select 2 matrix materials; material 1 : Boron, read from database, 
material 2. Aluminum, user supplied properties) 


♦MATRIX 

NMATX=2 

NM=1 MA=6 MAT=A IDB=y 

NM=2 MA=1 MAT=A IDB=n EL =E, v, G, a & 

Vl=a, b ,c ,d, ... K=c,j, k t % 


nr Note: The *a> k t data is only required if the ♦COND keyword has been 
used; where k a , k t are the axial and transverse thermal con- 
ductivities:, respectively. 
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is- Note: Required Format for User Supplied Isothermal Material Constants: 


Bodner-Partom: nemd = 1 

Elastic : Inelagtig: 

EL=E|_, Ej v A , Vy, G A , a A , a T Vl= D 0 , Zq, Zf, m, n,q 


Robinson Creep: nemd = 2 

Elastic : Inelastic : 

EL=E|_, Ej , Vy, G A , a A , ccy Vl= ctq, Eq, n, P, m, R, H 


Robinson Viscoplastic: ncmd= 3 
Elastic : Inelastic: 

EL=E l , E t v a , Vy, G a , a A , <Xy Vl= n, m, [L, k t , p, R, H, G 0 , q, to 


GVIPS: nemd =4 

Elastic : Inelastic: 

EL=E, v Vl= \l, k, R a , R k , B 0 , B v L 0 , L x , m, n, p, q, 

w, z 0 


Elastic Model: ncmd = 6 
EL=E L , E t v A , vy, G a , a A ,a T 

Table II. Non-isothermal Material Constants 

The Following materials are available for a non-isothermal analysis: 

MAC Materials Database 


Model 

Material 

Units 

mat 

Bodner-Partom 

Graphite 

MPa. sec, 0 C 

A 

nemd = 1 

Aluminum 

If u 

B 

GVIPS 

TIMETAL21S 

ksi, sec 

A 

nemd =4 





Required Format for User Supplied Non-isothermal Material Constants: 
Note: each of the following data statements are separate lines. 
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NTP=/ifpte 

TEM=7y, T 2 , T ntpts 
EAsE at 1 ’ E at 2 > — E at M pU 
ET=E TTi ’ E tt 2 ’ '■■ E TT n , pti 

NUAsV^jy ...V A J 

NUT^Vj’y. ) Vy*jt $ *“ 

GA=G^ r , Gat 2 ’ '"G 
ALPA=a AT , &at 2 ’" a AT 


nip Is 
ntpts 


nip is 


ALPTscxyy,, o. Tj , 


DO=D 0r , £?o r * 


nlpts 


zo=z 0r 

f> 

0 

N 

N 

0 

T \ 

1 4 2 


Z1=Z, 

lr l 

> z i T » 
1 r 2 

- Z l r 

BM=m 7 

V*r^ 

...m r 

* ntpts 

11 

z 

< 

, n T , . 

. . /It 

/ n/p/J 


<77y •••0r 1 


niplj 


® Note: For MF=99 (user defined material model) the following special for- 
mat is required 

*FIBER 

NF\BS= nfibs 


The following line is to be repeated for each fiber (nfibs): 

NF =nf f MS =ms f MF=99 f NPE=qpe EL=e,,e 2 ,...e npe & 
NP V=npv V\=v 1 ,v 2 ,...v npv K=k a , k t % 


^MATRIX 

NMATX=nmatx 


The following line is to be repeated for each matrix (nmatx): 

HM=nm m MS =ms m MA=99 m NPE -npe EL=e 7 ,e 2 .- -®npe & 
NPV=npv VI ^v^v& . Vnpv K=k a , k t % 
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where: 

ripe: 

total number of elastic constants 


a 7 1 ©2i • • 

elastic constants 
npv. 

total number of inelastic constants 


inelastic constants 
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4.2.12 RVE Data : 

Purpose: Select RVE representing desired fiber packing arrangement/architecture 

*MRVE 

\DP=idp 

«a* Note: The following data is entered on the same line as idp (except where 
noted) 

Without Interface : 

• For IDP = 1,2, or 3 

V F=vf 

• For IDP = 4 

VF=vf XA=xa 

• For IDP = 6 

VF =vf R =R 

• For IDP = 9 

VF1=vff RADI =rad1 VF2=vtf2 RAD2 =rad2 R =fl 

• For IDP = 11 

VF =vf1 RAD =rad1 R=f? 

With Interface: 

• For IDP = 1, 2,or 3 

VF=vf RAD=racf7 CPER -cper 

• For IDP = 4 

» Currently NOT Available 

• For IDP = 6 

® Currently NOT Available 

• For IDP = 9 

VF1 =vf1 RADI =radf CPER1 =cper1 MF2=vf2 & 

RAD2 =rad2 CPER2 =cper2 R=R 

• IDP = 11 

VF =vf1 RAD =rad1 R =R CPER=cper 


where: 

vf, vfl, vf2 = the fiber volume ratios 
rad, radl, rad2 = fiber radii 

cper, cperl, cper2 = ratios of interface thickness to fiber radius 
xa = length of the cross, see Fig. 6. 


46 



MAC M icromechanics A nalysis C ode 


= X/Y which defines the ratio of distances between fibers within a “ply” and 
those between a “ply” (see figure 7) 

idp:- ID of RVE as given below 

■»’ Note : For Double Periodic RVE: 

idp Description 

• 1 = Square Fiber, Square Pack (original 4-cell model) 

^ 1/(1 + A) 2 

RVE shown in Fig. 5 

• 2 = Square Fiber, Triangular (hexagonal) Pack 

V f £ 0.86602/(1+ A) 2 
RVE shown in Fig. 5 

• 3 = Square Fiber, Square Diagonal Pack 

Vy<;0.5/(l+A) 2 

RVE shown in Fig. 5 

• 4 = Cross Shaped Fiber, Square Pack 

V j <, 1 - 4(jc a) 2 

RVE shown in Fig. 5 

• 6 = Circular Fiber Approximation Rectangular or Square Pack 


v.s 08125 2 

if 

R> 1.0 

1 R(l+A) 



*(0.8125) 

if 

R< 1.0 

(1 + A) 2 




RVE shown in Fig. 5 
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• 9 = 


• 11 = 


• 99 = 


Two Different Size Square Fibers, Rectangular or Square 
Pack 




R i 


(1 + A 2 ) + (1 + Ai)- 


2,fRV t 




(l+A 2 ) + 


! /J 


(l+A,) 


RVE shown in Fig. 5 
Square Fiber, Rectangular Pack 

V f <. — = if R> 1.0 

1 R( l+A) 

V f <. if /?<1.0 

(l + A) 2 

RVE shown in Fig. 5 
User Defined RVE 

Example of RVE representing random packing shown in Fig. 
8. Required input shown in Appendix F. 


Note : For Triple Periodic RVE: 

• 99 = User Defined RVE 


«*■ Note: R, R f , R f2 , V fJ A = cper, A l = cperl and A 2 = cper2 are defined later in 
this subsection 

«*■ Note: The images of the RVE’s shown in Figs.5 and 8 were generated using 
the PATRAN option (sec 4.2.8) available in MAC. 


Example: Triangular packing with interface thickness 1% of fiber radius. 
3HVIRVE 

IDP=2 VF=35 RAD=0.07 CPER=0.01 % 
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IDP = 1 



IDP = 2 



IDP = 3 



IDP = 4 



IDP = 11 

Figure 5: RVE’s Available in MAC 
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Fig. 6 Cross Shaped Fiber Distance xa 



Fig.7 Hybrid Composite RVE; IDP = 9. Fiber Spacing Ratio, R=X/Y 
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er Note: Format for the user defined RVE’s, IDP = 99, see Fig. 8 for example. 


1) For 2J) RVE: (Each line of data must be on a separate line) 


NB=nb NG=ngr 

H —h-j, i •••• ftfig 

h' Inb 

CM=SS n b' ijSSnb'Z-SSnb.ng 


repeat data for the (nb x ng) 2-D RVE 


CM=SSi i ,SSi 2, SSf 


ng 


2) For 3 RVE: (Each line of data must be on a separate line) 


NA=na NB=nb NG=ng 


D = c/y , d 2 , ■ ■ ■ d na 
H=/Jji hg, — b n g 

h’ Inb 


CM=SS na i t i, SS na 2' i,--- SS na nb„ 1 

CM =SS na .i'i t i,SS na .i2,i’— ss na-l,nb„l 

repeat data for the (na x nbxngf) 3-D RVE 

CM =SS 1>1n g,SS 1i 2,ng<--- ss 1,nb,ng 


Where 

ss jj and ss/y ^ are the identifying material labels which are given in the 
following format: 

first character: F - for a fiber or M - for matrix 

second character: 1 , 2, 3, ... for fiber/matrix number 1 , 2, 3 ... 

na - number of subcells in the alpha direction 

nb - number of subcells in the beta direction 

ng - number of subcells in the gama direction 

h - the height of each subcell 

I - the length of each subcell 

d - the depth of each subcell 
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Example: (2-D user defined RVE (as shown above) 

*MRVE 

IDP=99 
NB=3 NG=3 

L-1.,1.,1. 

CMsM1.M1.M1 
CM=M1,M1,F1 
CM=F1,M1,M1 % 


is- Note: the fiber id. is FI and the matrix id. is Ml 
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Example: A 3-D User defined RVE (as shown above), representing a [0/90] 
continuous reinforced laminate composite. 

ra* Note: the fiber is material number 1 and the matrix is material number 2. 
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cross-section 1 cross-section 2 


*MRVE 

IDP=99 ASC=1 . 

NA=4 NB=2 NG=2 

L=1.,1. 

CM=M1,M1 rows 1-4 are for cross-section 1 

CM=M1,M1 

CM=M1,M1 

CM=F1,M1 

CM=M1 ,M1 rows 5-8 are for cross-section 2 

CM=F1,F1 

CM=M1,M1 

CM=F1 ,M1 % 
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4.2.13 Interface Data : (Optional, for IDP=1 ,2,3,9 or 11) 
Purpose: Specify interface layer properties 

INTERFACE 

NINT=n/nf 

The following line is to be repeated for each interface(nint): 
Nl =nm 1 M \=ncmd 1 MAJ=mat 1 IDB =idb 1 


% 

* 5 * Note: If idb=n see section 4.2.11 for format of additional input required. 
Where 

nint 

number of different interfaces 

nm: 

interface number 

ncmd: 

material model identifier for the interface 

mat 

material id letter 

idb\ 

= y use database properties 
= n user specified properties 

os- Note: Please refer to section 4.2.11 for ncmd and mat values 


Example: (Assuming one interface with user supplied material constants) 

•INTERFACE 

NINT=1 

Nl=1 Ml=4 MAT=C IDB=n & 

EL=1 1 700., 1 1 700..0.365, 0.365,4287.5,1 .,1 . & 

VI=0.8E-8,0.1 ,0.1 E-5,0.,0.85E-3, 0.05,1 .,1 .,1 .,3.3,1 .8,1 .35,1 .,0.01 % 
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4.2.14 Debond Data: (Optional) 

Purpose: Specify the subceil faces where debond can take place (see section 
3.6) 


2-D GMC Model Format 

*DEBOND 

Nll=n// 

NBI =nbi NGI =ngi FACE=nfc RN=rn BDN= bdstrn RS=rs & 
BDS=bdstrs T\=tmd 


% 

3-D GMC Model Format 

*DEBOND 

Nll=n// 

NAI=na/ NBI=nb/ NGI=ngi FACE =nfc RN =rn BDN=bdstrn RS=rs & 
BDS=bdstrs Tl =tmd 


% 

Where 

nii: 

number of subcell interfaces with debond 

nai, nbi and ngi: 

subcell indices (a, (3, y) 

nfc: 

used specify which “face” of the subcell is to be debonded 

For 2-D GMC 
= 1 Top Face 
= 2 Right Face 
= 3 Both Faces 

For 3-D GMC 
= 1 Top Face 
= 2 Right Face 
= 3 Back Face 


rn: 

R n value 
bdstrrr. 

debond stress, normal component (see Fig. 9) 


= 4 Top & Right Face = 7 ALL Faces 
- 5 Top & Back Face 
= 6 Right & Back Face 


rs: 

R s value 
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bdstrs : 

debond stress, shear component (see Fig. 9) 
tmd 

time at which debond criteria becomes active 
(defaults to tmd= 0. if not specified) 


Example 1 : 

debond subcell 1,1 top lace normal component 

Nll=1 

NBI=1 NGI=1 FACE=1 RN=1. BDN=100. 

Example 2: 

debond subcell 1,1 right face shear component 

Nl 1=1 

NBI=1 NGI=1 FACE =2 RS=1. BDN=100. 




FIBER 

(1,1) 



Top Face 


( 1 , 1 ) 


Right Face 


Figure 9: Subcell Faces For 2-D GMC 
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Figure 10: Subcell Faces For 3-D QMC 
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4.2.15 Plot Point Information : 

Purpose: Specify the frequency at which data will be written to output files for 
both X-Y data 

#CURVE 

NP=npmax 


% 

Where 

npmax: 

plot point increment 

Example:: (print out every fifth data point) 

*CURVE 

NP=5 % 
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4.2.16 Curve Data : 

a- Note: Currently a maximum of 5 curves may be specified per problem 

Specifying Curve Data For Macro (composite) Qu antities: 

*MACRO 

NT =nucuv 

Repeat the following line nucuv times: 

NC =nocu X=maidx Y=maidy NAM=tname 


% 

Where 

nocu: - curve number 


maidx and maidy variable options are: 
1 _ e n 7 - CT n 13 ' e n 

19-4 

22 - Total Time 

<N 

<N 

I 

CM 

8 " <^22 

14 ■ ^22 

20 - £22 

23 - Creep Time 

3 - 6-33 

9 “ 

15-233 

21 - 4 

24 - Temperature 

4 ‘ Y 23 

10- t 23 

16- Y 23 



5- Yis 

1 1 ' t 13 

1 7 ‘ Yl 3 



®’ Yl2 

12- t !2 

18 - Y'l2 




tname: name of plot file 


■or Note: The file(s) generated will be of the form fna/ne_mac.data. If the user 

desires to use the same tname for files 2-5, a double quote, ", is entered 
for tname. Those files then will have the form: 

tname 2_mac.data 
tname 3_mac.data 
etc. 

(see the example for more details) 
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EXAMPLE: 

♦MACRO 


NT=3 




NC=1 

X=1 

Y=7 

NAM=stress 

NC=2 

X=22 

Y=13 

NAM =“ 

NC=3 

X=23 

Y=13 

NAM= plot 


file created: stress_mac.data 
file created: stress2_mac.data 
% file created: plot_mac.data 


Specifying Curve Date for Micro (subcelh Quantities : 

♦MICRO 

UT=nucuv 

Repeat the following line nucuv times: 

NC =nocu CELL =nssel X=maidx Y=maidy NAM=tname2 

% 

Where 

nocu : - curve number 
nssel: - cell number 

tname 2 . - name of plot file 
miidx and miidy variable options are: 


1 - e n 

7 -an 

13-e'n 

19 -O u 

25 - «P U 

31 - Total Time 

2 _ e 22 

8 *a 22 

14-4 

20-022 

26-^22 

32 - Creep Time 

03 

■ 

w 

9 - a 33 

1 5 - e 33 

21 - <t > 33 

27 - * 33 

33 - Temperature 

4 - Y23 

1 0 - t 23 

16- Y23 

22 - 0 23 

28 - ^23 


5- Yi3 

11 -Tis 

17- y! 3 

23 - <t > 13 

29 - >P 13 


6 - Yi 2 

12-^12 

1 8 * y'i 2 

24 - 0 12 

30 - ^ 12 



«*■ Note: The quantities <t> and 'F are the possible internal state variables (consti- 
tutive model dependent). 
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«s* Note: The file(s) generated will be of the form fname2_mic.data. If the user 
desires to use the same tnameZ for files 2-5, a double quote, “, is 
entered for tname 2. Those files then will have the form: 

fname22_mic.data 
tname 23_mic.data 
etc. 

(see the example for more details) 

<>? Note: Subcell numbering is assigned according to the following algo- 
rithm. 

2- D case 

Do* IB=1,NB 
Do* IG=1,NG 

subcell number = NG*(IB-1 )+IG 

* continue 

3- D case 

Do * IA=1,NA 
Do* IB=1,NB 
Do * IG=1 ,NG 

subcell number = IA+(IB-1 )*NA+(IG-1 )*NA*NB 

* continue 

ca* Note: Refer to numbers in upper left corners of subcells in Figures on 
page 52 and 54 for examples of subcell numbering scheme. 


EXAMPLE: 

*MICRO 

NC=1 CELL=1 X=1 Y=7 NAM=cell % file created: cell_mic.data 
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5.0 Conclusion/Future Modifications 

A computationally efficient, user-friendly, comprehensive, micromechanics anal- 
ysis tool, MAC, has been presented that admits physically based viscoplastic 
deformation and life models, can analyze multiphased materials of interest in 
advanced propulsion systems, and can assist both the material scientist and 
structural analyst in developing, designing and analyzing strategic materials. How- 
ever, the development of this tool is far from complete. A number of future 
enhancements that are planned and currently underway include: 

• Provide an implicit integration algorithm to improve computational effi- 
ciency. 

• Provide the ability to input a multiaxial state of loading, be it stress or 
strain. 

• Provide the ability to input multiple runs within the same input data file. 

•The incorporation of damage evolution laws and failure criteria so as to 
provide life estimates. 


5.1 Acknowledgment 

The authors would like to thank Mr. Robert Goldberg for his diligent efforts and 
assistance in checking the accuracy and consistency of this manual by running 
numerous test cases. Future users will be, as the authors are, grateful. 
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6.0 Appendix A 

Sample Input File For A Mechanical Load Problem 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 


Load Type: 

Mechanical 

Load Component: 

22-direction (transverse to fiber) 

Load History: 

Cyclic 

Load Control: 

Strain 

Load History Data: 

e = 8.333X10 -4 , 


£ max = 0015 ■ 


*min = °- 


^initial = 0.0000024 

Micromechanics model: 

Double Periodicity 

Fiber Packing Arrangement: 

Hexagonal Pack at 35% fiber volume ratio 

Integration Algorithm: 

Predictor/Corrector 

Constituent Material Model: 

GVIPS 

Constituents: 

Fiber: SCS-6 

Matrix: TIMETAL21S 

Interface: fictitious weak interface for 

TIMETAL21S 

car Note: This example will take an extremely long time to run due to the elas- 
tic-perfectly viscoplastic definition of the fictitious weak Interface, 
as the overall time step is limited to that of the allowable interface 
time step. This slow computational response is NOT indicitive of 
GMC. 
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FOR MECHANICAL LOAD RUN 
•PRINT 
NPL=0 % 

•LOAD 

LCON=2 LOP=2 LSS=1 % 

*MECH 

NPTW=2 Tl=0.,18. L0=0.,0.015 % 

*TREF 

TREF=650. % 

•MODEL 
MOD=1 % 

•SOLVER 

NTF=2 ISTM=0.0000024 ERR=0.1E-2 % 

•FIBER 

NFIBS=1 

NF=1 MF=4 MAT=B IDB=y % 

•MATRIX 

NMATX=1 

NM=1 NM=4 MAT=A IDB=y % 

*MRVE 

IDP=2 VF=0.35 RAD =0.07 CPER=0.1 % 

•INTERFACE 

NINT=1 

Nl=1 Ml=4 MAT=C IDB=n & 

EL=11700.,11700.,0.365, 0.365,4287.5,1., 1. & 

VI=0.8E-8,0.1 ,0.1 E-5.0..0.85E-3, 0.05,1 .,1 .,1 .,3.3,1 .8,1 .35,1 .,0.01 % 
•CURVE 
NP=2 % 

•MACRO 

NT=1 

NC=1 X=2 Y=8 NAM=apdxa % 

•END 
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The following figure was obtained from the x-y plot data file produced by the 
present example. 



e 


XX 


(%) 


«*■ Note: It is recommended that a new user construct a mac input file using the 
data given in this appendix and then check to see if the same result plot 
is obtained. 
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7.0 Appendix B 

Sample Input File For A Thermal L oad Problem 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 

Thermal 

T = 0 . 01 , 


Load Type: 

Load History Data: 


T = 371 1 

1 max 


T - = 21.1 

x min 


^initial = 17.505 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 
Constituents: 


Double Periodicity 

Square Pack at 35% fiber volume ratio 

Predictor/Corrector 

Bodner-Partom 

Fiber: Graphite 
Matrix: Aluminum 


cr Note: This problem is taken from reference 1 pg.238 
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FOR THERMAL LOAD RUN 
•PRINT 
NPL=0 % 

•LOAD 
LCON=1 % 

THERM 

NPTT=3 Tl=0.,35010.,70020. & 

TE=371. 1,21. 1,371.1 % 

TREF 

TREF=371.1 % 

•MODEL 
MOD=1 % 

•SOLVER 

NTF=2 ISTT=1 7.505 ERR=10.E-2 % 
•FIBER 
NFIBS=1 

NF=1 MF=1 NDPT=2 MAT=A IDB=y % 
•MATRIX 
NMATX=1 

NM=1 MM=1 NDPT=2 MAT=B IDB=y % 
•MRVE 

IDP=1 VF=0.30 % 

•CURVE 
NP=5 % 

•MACRO 

NT=1 

NC=1 X=24 Y=1 NAM=thermal % 

•END 
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The following figure was obtained from the x-y plot data file produced by the 
present example. 


xlO ' 2 



S3- Note: It is recommended that a new user construct a mac input file using the 
data given in this appendix and then check to see if the same result plot 
is obtained. 
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8.0 Appendix C 

Sample Input File For A Laminate Problem 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 

Strain Control 

e = 0.01 , 


Load Type: 

Load History Data: 


= ooi , 


A'cmut.' = 0.00025 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 
Constituents: 


Laminate with Double Periodicity 

Square Pack at 46% fiber volume ratio 

Forward Euler 

Bodner-Partom 

Fiber: Boron 
Matrix: Aluminum 
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test of laminate strain control 
•PRINT 

NPL=0 % 

•LOAD 

LCON=2 LOP=1 LSS=1 % 

•MECH 

NPTW=2 Tl=0.,1. L0=0.,0.01 % 

•MODEL 

MOD =3 MATSYS=1 NLY=1 THK=1. ANG=0. % 
•SOLVER 

NTF=1 NPTS=2 TIM=0.,1 . STP=0.00025 % 
•FIBER 

NFIBS=1 

NF=1 MS=1 MF=6 MAT=A IDB=y % 

•MATRIX 

NMATX=1 

NM=1 MS=1 MM=1 MAT=A IDB=y % 

*MRVE 

IDP=1 VF=0.46 % 

•CURVE 

NP=10 % 

•MACRO 

NT=1 

NC=1 X=1 Y=7 NAM=apdxc % 

•END 
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The following figures were obtained from the x-y plot data file produced by the 
present example. 



From reference 1 Fig. 8-5, pg. 235 



From reference 1 Fig. 8-7, pg. 237 
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9.0 Appendix D 


Sample Input File For A Laminate Problem 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 

Load Type: Stress Control 

Load History Data: d = 0.01 ksi/sec, 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 
Constituents: 


^ constant 1-125 

[±45] s Laminate with Double Periodicity 

Square Pack at 45% fiber volume ratio 

Forward Euler 

Bodner-Partom 

Fiber: Boron 
Matrix: Aluminum 


«*■ Note: This problem is taken from the reference: 

Analysis of Metal-Matrix Composite Structures-ll. 

Laminate Analysis, Arenburg, R. T. and Reddy, J. N., 
Computers and Structures, Vol. 40, N.6, pp. 1369-1385, 1991 . 


74 


MAC: Micromechanics A natysis C ode 


test of [+-45]2s laminate strain control 
•PRINT 
NPL=0 % 

•LOAD 

LCON=2 LOP=1 LSS=2 % 

*MECH 

NPTW=2 Tl=0.,4500. LO=0.,45. % 

TREF 
TREF=0. % 

•MODEL 

MOD=3 MATSYS=1 NLY=4 THK=0.25, 0.25, 0.25,0.25 ANG=45.,-45.,-45.,45. % 
•SOLVER 

NTF=1 NPTS=2 TIM=0.,4500. STP=1.125% 

•FIBER 

NFIBS=1 

NF=1 MS=1 MF=6 MAT=A IDB=n & 

EL=58.E3, 58.E3, 0.20, 0.20, 24.1 7E3,6.3E-6,8.28E-6 % 

•MATRIX 

NMATX=1 

NM=1 MS=1 MM=1 MAT=A IDB=n & 

EL=9.53E3,9.53E3, 0.33, 0.33, 3.58E3, 21 .06E-6.21 .06E-6 & 
VI=1.E4,49.,63.,300.,4.,1. % 

•MRVE 

IDP=1 VF=0.45 % 

•CURVE 
NP=10 % 

•MACRO 

NT=1 

NC=1 X=1 Y=7 NAM=apdxd % 

•END 
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The following figure was obtained from the x-y plot data file produced by the 
present example. 



From Arenburg and Reddy, Fig. 16, pg. 1382. 
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10.0 Appendix E 


Sample Input File For Triple Periodic GMC 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 


Load Type: 

Load History Data: 


Strain Control 

e = 0.01 , 

^max = 001 • 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 
Constituents: 


^constant = 0.00025 

Triple Periodic GMC 

Square Pack at 46% fiber volume ratio 

Forward Euler 

Bodner-Partom 

Fiber: Boron 
Matrix: Aluminum 


This problem uses the 3-D GMC model to simulate a [0/90] laminate as shown 
below. 
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test of gmc3d model 0/90 laminate simulation 
•PRINT 
NPL=0 % 

•LOAD 

LCON=2 LOP=2 LSS=1 % 

*MECH 

NPTW=2 Tl=0.,1. L0=0.,0.01 % 

•MODEL 
MOD=2 % 

•SOLVER 

NTF=1 NPTS=2 TIM=0.,1. STP=0.00025% 
•FIBER 
NFIBS=1 

NF=1 MF=6 MAT =A IDB=y % 

•MATRIX 

NMATX=1 

NM=1 MM=1 MAT=A IDB=y % 

*MRVE 

IDP=99 ASP=1 . 

NA=4 NB=2 NG=2 

D=0.67823, 0.321 77,0.67823,0.321 77 

H=0.67823, 0.321 77 

L=0.321 77,0.67823 

CM=2,2 

CM=2,2 

CM=2,2 

CM=1,2 

CM =2,2 

CM=1,1 

CM=2,2 

CM=1 ,2 % 

•CURVE 
NP=10 % 

•MACRO 

NT=1 

NC=1 X=2 Y=8 NAM=apdxe % 

•END 
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The following figure was obtained from the x-y plot data file produced by the 
present example. 



■y Note: It is recommended that a new user construct a mac input file using the 
data given in this appendix and then check to see if the same result plot 
is obtained. 
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11.0 Appendix F 

Sample Input File For A Us er Defined RVE 
The following example is used to explain the control blocks in more detail. 
Problem Summary : 

Load Type: Strain Control 

Load History Data: e = 0.01 , 

z max = 0001 . 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 
Constituents: 


^constant 0-1 

Double periodic 

Random composite (see figure 8) 

Forward Euler 

Bodner-Partom 

Fiber: Boron 
Matrix: Aluminum 
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random composite rve using 2-d gmc 
•PRINT 
NPL=0 % 

•LOAD 

LCON=2 LOP=1 LSS=2 % 

•MECH 

NPTW=2 Tl=0.,0.1 LO=0., 0.001 % 

•MODEL 
MOD=1 % 

•SOLVER 

NTF=1 NPTS=2 TIM=0.,0.1 STP=0.1 % 

•FIBER 

NFIBS=1 

NF=1 MF=6 MAT=A IDB=y % 

•MATRIX 

NMATX=1 

NM=1 MM=1 MAT=A IDB=y % 

*MRVE 
IDP=99 
NB=14 NG=14 

CM=F1,M1,M1,M1 , Ml ,M1, FI, FI, Ml, Ml, FI.FI.MI.FI 
CM=F1 , Ml, FI, F1.M1.M1, FI .FI.MI.MI, Ml .M1.M1.F1 
CM=M1 .M1.F1.F1, Ml, Ml, Ml, Ml, Ml, M1.F1.F1, Ml, Ml 
CM*M1 , Ml , Ml ,M1, MI.MI.MI.FI, FI.MI.FI.FI, Ml, Ml 
CM=M1,M1,F1,F1,M1,M1,M1,F1,F1,M1,M1,M1,M1,M1 
CM*F1, Ml, FI, FI, Ml, Ml, Ml, Ml, Ml, Ml, F1.F1.M1 .fi 
CM=F 1 ,M1 ,M1 ,M1 ,M1,F1 ,F1 ,M1 ,M1 ,M1 ,F1 ,F1 ,M1 ,F1 
CM=M1 ,M1,M1 ,M1 ,M1 ,F1 ,F1 ,M1 ,M1 ,M1 ,M1 ,M1 ,M1 ,M1 
CM=M1, FI, FI, Ml, Ml, Ml. Ml, Ml, FI, FI, Ml, FI, FI, Ml 
CM=M1, FI, FI, Ml, FI, FI, MI.MI.FI.FI, Ml, FI, FI ,M1 
CM=M1, Ml, Ml, Ml, FI.FI.FI.FI, Ml, Ml, Ml, Ml, Ml, Ml 
CM=F1, FI, Ml, Ml, Ml, Ml, FI, FI, Ml, Ml, Ml, FI, FI, Ml 

CM=F1, FI, Ml, Ml, Ml, Ml .Ml, Ml, Ml, Ml, Ml, FI, FI, Ml 

CM=M1, Ml, Ml, Ml, Ml, Ml, Ml, Ml, Ml, Ml, FI, FI, Ml, Ml % 

•CURVE 
NP=1 % 

•MACRO 

NT=1 

NC=1 X=1 Y=7 NAM=apdxf % 

•END 


81 



MAC Mkxomedharics A nalysis Code 



Fiber Packing Arrangement: Square Pack at 35% fiber volume ratio 

Integration Algorithm: Forward Euler 

Constituent Material Model: Elastic 
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test of biaxial load 
•PRINT 
NPL=7 % 

•LOAD 

LCON=2 LOP =7 LSS=1 % 

•MECH 

NPTW=2 Tl=0.,1.5 L0=0.,0.015 
NPTW=2 Tl=0.,1 .5 LO=0., 0.010 % 

TREF 

TREF=900. % 

•MODEL 
MOD=1 % 

•SOLVER 

NTF=1 NPTS=2 TIM=0.,1 .5 STP=0.015 % 

•FIBER 

NFIBS=1 

NF=1 MF=4 MAT=B IDB=n & 

EL=58.E3, 58. E3, 0.32,0.32, 23.2E3.2. 1 E-06,2. 1 E-06 & 

VI=0.8E-9,1 .E20, 0.1 E-5,0.,0.85E-4, 0.05,1 .,1 .,1 .,3.3,1 .8,1 .35,1 .,0.01 % 
•MATRIX 
NMATX=1 

NM=1 MM=4 MAT=A IDB=y % 

•MRVE 

IDP=1 VF=0.35 % 

•CURVE 
NP=1 % 

•MACRO 

NT=2 

NC=1 X=1 Y=7 NAM=biaxl 
NC=2 X=2 Y=8 NAM=biaxt % 

•END 
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13.0 Appendix H 

Sample Input File For A User Defined Material Model 


The following example is used to explain the control blocks in more detail. 
Problem Summary : 

Load Type: Stress Control 

Load History Data: d = 0.01 ksi/sec, 


Micromechanics model: 
Fiber Packing Arrangement: 
Integration Algorithm: 
Constituent Material Model: 


^constant ~ 1-125 

[±45] s Laminate with Double Periodicity 
Square Pack at 45% fiber volume ratio 
Forward Euler 

Bodner-Partom using user subroutine 
USRMAT 


Constituents: 


Fiber: Boron 
Matrix: Aluminum 


*r Note: This problem is the same as that given in Appendix D. The purpose here 
is to demonstrate how to use the USRMAT option. 
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use of usrmat for a [+-45]2s laminate strain control problem 
'PRINT 
NPL=0 % 

♦LOAD 

LCON=2 LOP=1 LSS=2 % 

*MECH 

NPTW=2 Tl=0.,4500. LO=0.,45. % 

♦TREF 
TREF=0. % 

♦MODEL 

MOD=3 MATSYS=1 NLY=4 THK-0.25,0.25,0.25,0.25 ANG=45.,-45.,-45.,45. % 
♦SOLVER 

NTF=1 NPTS=2 TIM=0.,4500. STP=1.125 % 

♦FIBER 

NFIBS=1 

NF=1 MS=1 MF=99 NPE=7 EL-SS.ES.SS.ES.O^.O^^.-^ES.e.SE-e.S^SE-e & 
NPV-6 VI-0.,1.,1.,1.,1.,1. % 

♦MATRIX 

NMATX-1 

NM=1 MS=1 MM =99 NPE=7 & 

EL=9.53E3 I 9.53E3,0.33,0.33 I 3.58E3,21 .06E-6.21 .06E-6 & 

NPV=6 Vl=1 .E4 I 49.,63.,300.,4.,1 . % 

♦MRVE 

IDP=1 VF=0.45 % 

♦CURVE 
NP=10 % 

♦MACRO 

NT=1 

NC=1 X=1 Y=7 NAM=apdxg % 

♦END 
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The USRMAT subroutine 


The USRMAT subroutine is used here to implement the Bodner-Partom Visco- 
plastic model currently available in MAC 



SUBROUTINE USRMAT(DSA,S A,PEJ>V, 

& TIME,TSTEP,CTEMP,DTEMPR,TREFJ4IO,NPE,NPV,NCE) 


c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


purpose: user material model 

note: 1) in this subroutine, [sa] and [dsa] contain the 

“micro” quantities for aboudi’s micromechanics model 

2) arrangement of [dsa] & [sa] arrays 


variable location 


1 

1 strain rate 
1 

(1-6) 

1 stress rate 
1 

(7-12) 

1 inelastic 

1 strain rate 
1 

(13-18) 

1 12 “slots” 

(19-30) 

1 for state variables 


I thermal strain rate (31-36) 


IMPLICIT DOUBLE PRECISION (A-H,0-Z) 
DIMENSION DSA(36),SA(36),PE(NPE),PV(NPV) 


c 

c default user supplied material routine 



c 



DIMENSION SS(6),S(6),R(6) 


c 
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c write(nio,*) ‘bodner: ctemp ‘,ctemp,’ time ‘,time 

c 

c extract appropriate viscoplastic material constants 

DO = PV(1) 

ZO = PV(2) 

Z1 = PV(3) 

BM = PV(4) 

AN = PV(5) 

Q = PV(6) 

c 

c copy stress from [sa] to [s] 

c 

SS(1) = SA(7) 

SS(2) = SA(8) 

SS(3) = SA(9) 

SS(4) = SA(10) 

SS(5) = SA(ll) 

SS(6) = SA(12) 

c 

c compute the deviatoric stress [s] in the subcell 

c 

TEMP = (SS(1) + SS(2) + SS(3))/3. 

S(l) = SS(1) - TEMP 
S(2) = SS(2) - TEMP 
S(3) = SS(3) - TEMP 
S(4) = SS(4) 

S(5) = SS(5) 

S(6) = SS(6) 


c 

AJ2=0.5*(S(1)**2+S(2)**2+S(3)**2)+S(4)**2+S(5)**2+S(6)**2 
SQ3AJ = DSQRT( SS(1)**2 + SS(2)**2 + SS(3)**2 + 

& 2*(SS(4)**2+SS(5)**2+SS(6)**2) ) 

SQ2=1. 414215 
EF(SQ3AJ.EQ.O.) THEN 
CALL ZEROR(R,6) 

ELSE 

R(l) = SS(1)/SQ3AJ 
R(2) = SS(2)/SQ3AJ 
R(3) = SS(3)/SQ3AJ 
R(4) = SQ2*SS(4)/SQ3AJ 
R(5) = SQ2*SS(5)/SQ3AJ 
R(6) = SQ2*SS(6)/SQ3AJ 
END IF 
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c if d0=0 then assume elastic and zero-out 
c [dsa( 13-30)] (inelastic strain rate and 
c internal variable rates), then return 

c 

IF(DO.EQ.O) THEN 
DO 100 JJ=13,30 
DS A(JJ) = 0.0 
100 CONTINUE 
RETURN 

c 

c inelastic 

c 

ELSE 


ZEF = ZO + Q*SA(20) + 

& ( 1 -Q)*(R( 1 )*S A(2 1 )+R(2)*S A(22)+R(3)*S A(23)+ 

& R(4)*SA(24)+R(5)*SA(25)+R(6)*SA(26)) 


C 

IF(AJ2 .EQ. 0.) THEN 

AL=0.0 

ELSE 

ARG1=ZEF**2.0/(3.*AJ2) 

CON= 5*(AN+ 1 .)/AN 
ARG=CON*(ARGl)**AN 
AL=D0/(DEXP(ARG)*DSQRT(AJ2)) 
ENDIF 


c inelastic strain rates 


DSA(13) = AL*S(1) 
DSA(14) = AL*S(2) 
DSA(15) = AL*S(3) 
DSA(16) = 2*AL*S(4) 
DSA(17) = 2*AL*S(5) 
DSA(18) = 2*AL*S(6) 


c plastic work rate 

C WPD = S(1)*DSA(13) + S(2)*DSA(14) + S(3)*DSA(15) + 
& S(4)*DSA(16) + S(5)*DSA(17) + S(6)*DSA(18) 



c state variable rates 
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DSA(19) = WPD 
Z0M=BM/Z0 
ZD=Z0M*(Z 1 -ZEF)*WPD 
DSA(20)=ZD 


DSA(21)=ZD*R(1) 

DSA(22)=ZD*R(2) 

DSA(23)=ZD*R(3) 

DSA(24)=ZD*R(4) 

DSA(25)=ZD*R(5) 

DSA(26)=ZD*R(6) 

ENDIF 

RETURN 

END 
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